\chapter{计量经济学理论与软件应用}
\section{回归模型}
\subsection{简单的回归分析}
若模型设置如下：
\begin{align*}
 lnwage&=\beta_0+\beta_1 educu \\
 lnwage&=\beta_0+\beta_1 educ+\beta_2 exper+\beta_3tenure+u 
\end{align*}
代码可如下撰写：
\begin{lstlisting}
library(wooldridge)
data("wage1")
lnwageModel <- lm(lwage ~ educ, data = wage1)
lnwageModel <- lm(lwage ~ educ+ exper + tenure, data = wage1)
summary(lnwageModel)
\end{lstlisting}

\subsection{公式的通俗写法}
\includegraphics[scale=0.7]{formula.png}

\subsection{假设检验：系数约束检验(\texttt{car})}
对于模型：
\[  prestige =\beta_0+\beta_1income + \beta_2education+u \]

若想检验假设$ H_0:\beta_1=\beta_2 \text{以及} H_0:2\beta_1=\beta_2$，则可如下撰写代码：
\begin{lstlisting}
library(car)
data(Duncan)
mod.duncan <- lm(prestige ~ income + education, data=Duncan)
linearHypothesis(mod.duncan, "income = education") # for  H0:beta1=beta2 
linearHypothesis(mod.duncan, "2*income = education") # for H0:2beta1=beta2 
\end{lstlisting}



\subsection{其他检验}
\begin{itemize}
	\item  自相关检验: 注意数据的顺序是排好的，因为自相关检验针对的时间序列。
	\begin{lstlisting}
library(lmtest)
dwtest(States$pop~States$pay+States$region)
	\end{lstlisting}	
	\item  异方差检验:
	\begin{lstlisting}
library(lmtest)
bptest(States$pop~States$pay+States$region)
\end{lstlisting}
	\item 多重共线性：
	\begin{lstlisting}
data("States") # `基础包数据集'
lmout <- lm(pop~pay+region,data=States)
vif(lmout) %>% sqrt()>2  # `膨胀因子'
	\end{lstlisting}
	
\end{itemize}


\subsection{时间序列：单位根检验}
library(fUnitRoots)
\begin{itemize}
	\item \texttt{adfTest}: Augmented Dickey--Fuller test for unit roots,
	\item \texttt{urersTest}: Elliott--Rothenberg--Stock test for unit roots,
	\item \texttt{urkpssTest}: KPSS unit root test for stationarity,
	\item \texttt{urppTest}: Phillips--Perron test for unit roots,
\item \texttt{urspTest}: Schmidt--Phillips test for unit roots,
\item \texttt{urzaTest}: Zivot--Andrews test for unit roots.
\end{itemize}
最常用的ADF检验命令：

\lstinline|adfTest(x, lags = 1, type = c("nc", "c", "ct"), title = NULL,description = NULL)|

\begin{lstlisting}
x = rnorm(1000)
# A time series which contains a unit-root:
y = cumsum(c(0, x))
## adfTest - 
adfTest(x)
adfTest(y)
\end{lstlisting}

做时间序列建模经常使用\lstinline|forcast|包。
\begin{itemize}
	\item \lstinline|checkresiduals()|就是检查白噪声的函数。LB检验的原假设是序列是独立的。
\end{itemize}

\subsection{虚拟变量回归：因子变量}
多类别的变量，作回归时，可直接输入，自动生成虚拟变量，不必担心虚拟变量陷阱。如，
\[ pop=\beta+\beta_1 pay+\beta_2 region+u \]

其中\texttt{region}就是虚拟变量。在确保\texttt{region}是因子变量前提下，可直接如下撰写代码：
\begin{lstlisting}
data("States") # `基础包数据集'
lmout <- lm(pop~pay+region,data=States)
summary(lmout)
\end{lstlisting}


\subsection{Probit 或者 Logit 模型}
\subsubsection{定义}
Logit模型: 
\[   P(y=1|X=x)=\frac{\exp(\bm{x}'\beta)}{1+\exp(\bm{x}'\beta)} \]

Probit模型:
\[ P(y=1|X=x)= \frac{1}{\sqrt{2\pi}}\exp(-\frac{\bm{x}'\beta}{2})\]

\subsubsection{估计}
Probit/Logit模型属广义线性模型中，使用\texttt{glm}命令：

\lstinline|glm(formula, family = gaussian, data, ...)|

评论：
\begin{itemize}
	\item 与普通线性回归一个最大的差别就是多了\texttt{family}参数
	\item \texttt{family}是对误差分布的设置，对于Probit/Logit而言，是二项式分布，使用前面关于婚外情的\hyperlink{2-1}{因子数据}：\\
\lstinline|glmout <- glm(ynaffair~gender+age+yearsmarried,data = Affairs,family = binomial(link = 'logit'))|
\end{itemize}


\subsubsection{Probit/Logit模型的解释}


\paragraph{拟合优度} 一般使用伪$ R^2$:
\[R^2=1-\frac{\mathscr{L}_{ur}}{\mathscr{L}_0}\]

	其中：$ \mathscr{L}_{ur} $是估计的对数似然值，$ \mathscr{L}_0 $是只有截距项的0模型的对数似然值。一般不报告，要使用\texttt{logLik(glmout)}抽取对数似然函数。
	
	\paragraph{系数的解释}：Probit/Logit模型可以写为$ P(y=1|\bm{x})=G(\bm{x\beta}) $，其中$G$为对数单位分布函数或标准正态分布函数，其偏效应为，
	\[ \frac{\Delta \hat P(y=1|\bm{x})}{\Delta x_j}=g(\bm{x\bm\beta})\hat \beta_j \]
	
	可见在不同的$ \bm{X}$上，偏效应是不同的。多数时候，取$ g(\bm{X\beta}) $的样本平均作为替代，此时，关键在于密度函数$ g(\bm{X\beta}) $的计算，一般使用\lstinline{predict}即可。
	
	要注意\lstinline|predict(object, newdata = NULL,type = c("link", "response", "terms"))|中\lstinline|type|的选择，缺省时给出的是对数odds比(\lstinline|"link"|)，\lstinline|"response"|才给出概率估计。

\begin{verbatim}
library(wooldridge)
data("mroz")
glmout <-  glm(inlf~nwifeinc+educ+exper+I(exper^2)+age+kidslt6+kidsge6,family = binomial(link = 'logit'),data = mroz)
summary(glmout)
predict(glmout,type = 'response')	
\end{verbatim}

\section{与stata的交互}
步骤如下：
\begin{itemize}
	\item 安装\texttt{RStata}包；
	\item 在R的启动环境文件\texttt{Rprofile.site}中设置\texttt{stata}的安装路径，即增加一行，如我的就是：
	\lstinline{options(RStata.StataPath = "\"D:\\Program Files (x86)\\Stata14\\StataMP-64\"")}；
\end{itemize}

在R中调用stata的一个例子：
\begin{figure}[H]
\begin{lstlisting}
library(RStata)
stata_src <- '
set more off
sysuse auto
reg mpg weight
'
stata(stata_src,stata.version=14) #`输命令'
stata('foo.do') #`直接调用do文件'
\end{lstlisting}	
\end{figure}